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1.0 


Introduction 


The major task carried out under Contract 
NAS2-5876 was the development of theoretical and 
computational techniques for calculating the time depen- 
dent electromagnetic response of a radially inhomogeneous 
moon. The importance of this task was two-fold. The 
techniques were first used to analyze the experimental 
data from the LSM (Lunar Surface Magnetometer) thus pro- 
viding the first in-depth diagnostics of the Lunar 
interior (Sonett et. al., 1971, a and b) . The theory 
was also incorporated into an existing computer code 
designed to calculate the thermal evolution of planetary 
bodies. The previous calculations (Sonett et. al., 1969) 
had relied only on unipolar induction. The new program 
will provide a tool for examining the effect of heating 
from the TE mode (poloidal magnetic field) as well as 
the TM mode (toroidal magnetic "field) . 

The theoretical discussion of the problem of 
the electromagnetic response of a radially inhomogeneous 
moon is presented in Section 2. The method of obtaining 
various potential solutions for the vector wave equation 
is shown along with the spherical harmonic expansion of 
these potentials. Formulas, are' derived for the radial 
equations which reduce the problem of obtaining the 
solutions to an initial value calculation. 
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The theory is applied to the problems of the 
near surface thermal gradient of the moon, and the 
overall lunar electrical conductivity in Sections 3 and 
4, respectively. By comparing the experimental data with 
calculated results for lunar models, the thermal gradient 
is bracketed between 1.5 and 4°K/km near the lunar sur- 
face. A computer inversion using the complete experimen- 
tal spectrum is invoked in Section 4 to obtain the lunar 
conductivity. This work has been reported in Sonett 
et. al. (1971, a and b) . The major discovery in this 
work is a prominent peak in the electrical conductivity 
at a radius of 1500 kilometers. The lunar conductivity 
exhibits a precipitous fall off on -either side of 1500 
kilometers. 

A new version of the computer code used to 
generate the data for Sonett et. al. (1969) is described 
in Section 5. The major change is the inclusion of the 
time dependent' electric fields for the joule heating. 

The program with all of its subroutines has been checked 
and is operational from an experimental aspect. Because 
of the inherent complexity of the overall heating problem, 
more short computer runs must be laid out and made in 
order to determine an optimum mode of investigation for 
this program. 
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2.0 


Electromagnetic Induction in the Moon 


2 . 1 Mechanism 

Electromagnetic induction -in the Moon is depen- 
dent upon the existence of a large class of discontinuities 
and waves in the solar wind. The forcing function which 
drives the induction is composed of .the steady interplane- 
tary magnetic field upon which is superimposed a hydromag- 
netic radiation continuum due to (1) plasma waves arising 
presumably in the solar atmosphere and subsequently 
convected and propagated outwards, (2) waves due to local 
instabilities in the solar wind, and (3) discrete large 
amplitude events such as collision-free shock waves, 
tangential discontinuities, and Alfven waves. 

The formal statement of the electromagnetic 

interaction problem rests upon division of the fields 

into parts which correspond to toroidal and poloidal 

magnetic fields. These are excited respectively by a) 

the interplanetary electric field (as perceived in a 

frame of reference comoving with the Moon) given by 

E = v x B , where v is the velocity of the Moon rela- 
— m — — — 

tive to the solar wind and _B the instantaneous inter- 
planetary magnetic field, and b) the time rate of change 
of the interplanetary field, B . The two modes corres- 
pond respectively to transverse magnetic (TM) and transverse 
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electric (TE) excitation. Both modes display strong 
frequency dependence (f = uj/2tt where oo is the angular 
frequency in radian/sec) . The magnitude of the TE 
transfer function goes to unity as frequency, £ -* 0 , 
and increases with increasing, f . The TE mode currents 
which close wholly in the lunar interior, tend to be 
concentrated where the waves damp substantially. With 
increasing frequency, the poloidal magnetic field becomes 
compressed into shells of decreasing thickness. The 
exclusion of the field from the interior of the moon, 
coupled with the confinement of the field at the lunar 
surface by the solar wind pressure on the sunlit side of 
the moon produces a magnification of the resultant 
magnetic field at the lunar surface. 

The TM mode attains peak value for f = 0 ; it 
remains approximately constant with increasing f until 
a combination of core and lithospheric constitutive param- 
eters forces the currents to pass wholly through the crust 
whereupon it decreases with further increase of f 
This mode is responsible for the steady state bow wave 
phenomena . 

Theoretical treatments of the lunar interaction 
with the solar wind generally include the effect of the 
solar wind dynamic pressure in confining the induced field 
lines (Blank and Sill, 1969; Schubert and Schwartz, 1969) . 



This is provided in the model by a field confining surface 
current layer in the solar wind just ahead of the lunar 
surface. Preliminary examination of the lunar response 

using the LSM data shows a strong amplification of inci- 

r 

dent tangential discontinuities whose free stream proper- 
ties are monitored by Explorer 35 (Sonett et. al . , 1971 a and b) 
The amplification occurs only for the vector components 
tangential to the surface; the normal component tends to 
follow the interplanetary value. Thus, the existence of 
a thin confining current layer appears verified for the 
sunward side of the Moon. The very strong excitation 
implies that the lines of force are confined within the 
Moon to a crustal layer having an electrical conductivity 
substantially less than that of the deeper layers. 

2 . 2 Theory of the Induction for Harmonic Excitation 

A complete treatment. of the electromagnetic 
interaction would require appropriate matching of the 
interior fields to those in the plasma surrounding the 
Moon. A less complicated interaction model can be used 
for the present day moon-solar wind electromagnetic induc- 
tion problem because the flow field exterior to the sunlit 
side of the moon shows no marked perturbation from the 
presence of the Moon aside from the diamagnetic cavity 
(Colburn et. al., 1967; Ness et. al., 1967). Based on 



this lack of a scattered field up stream in the plasma, 
an analysis for the sunward hemisphere is carried out for 
an inhomogeneous Moon immersed within a perfectly con- 
ducting space. The compressional effect of the solar wind 
momentum flux is taken into account mathematically by 
postulating an electrical current layer just above the 
lunar surface. 

Previous discussion of the electromagnetic 
response of the Moon to solar wind forcing functions have 
been limited to models of concentric, constant electrical 
conductivity spherical shells (Blank and Sill, 1969, 
Schwartz and Schubert 1969, Schubert and Schwartz, 1969, 
Sill and Blank, 1970) . Also, the analyses Of Blank and 
Sill (1969) . and Sill and Blank (1970) are valid only for 
low frequencies.. Their work corresponds to the retention 
of only the first term in the spherical harmonic expansion 
of the forcing field. Electromagnetic waves of frequency 
high enough to be confined to the near surface region of 
constant temperature gradient require for their complete 
description higher order terms in a spherical harmonic 
expansion. The theory presented in this section is appli- 
cable to a Moon whose electrical conductivity, a , 
varies with radial position in an arbitrary manner. In 
addition, higher order terms in a spherical harmonic 
expansion of the fields are retained throughout the 


- 6 - 


solar wind 


theoretical development. 

We consider the interaction of a 
forcing field of the form 


H = y H q exp <j i - cut) 


( 1 ) 


E = x |ivH Q exp i (-^p - cut) 


(2) 


with a radially inhomogeneous Moon. The cartesian 
coordinate system (x, y, z) with unit vectors x, _y, 
is fixed relative to the Moon and has its origin at the 
Moon's center. The Moon moves with speed v in the 
negative z-direction. The quantities H Q and A = v/f 
are the amplitude and wavelength, respectively, of the 
magnetic field oscillation and the circular frequency co 
is 2irv/A . The amplitude of the electric field' 
oscillation is pvH Q , where p is the permeability of 
the solar wind. We follow the formalism of Wyatt (1962) 
who considered the scattering of electromagnetic plane 
waves from inhomogeneous spherically symmetric objects 

The following theoretical development, as applied to the 
Moon problem, was first presented by Schubert and Schwartz 
(1970) . This work was then incorporated into the LSM 
report (Sonett et. al., 1971 a,b) . 
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in connection with the problem of determining the nuclear 
scattering cross-sections. His work appears to have been 
based on earlier work of Lahiri and Price (1939) . 

2.2.1 Formal Solution for a(r) 

The solutions of Maxwell's equations can be 
represented as the superposition of two linearly indepen- 
dent fields, the transverse electric (JS m , H™) and the 
transverse magnetic (_E , H ) . It can be shown that these 

two sets of fields are completely determined from two 
potentials, n and fi which satisfy the equations 

v 2 fi m + k 2 n m = 0 (3) 


^ £ < rnS > + k2[f - 0 


(4) 


where 


2 2 

k = (jo Pe + lapcu 


(5) 


k^ = icoe - <7 


(6) 
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and p and e are the magnetic permeability and elec- 
trical permittivity of the Moon. Although some possibility 
exists that local effects are significant in increasing 
the value of p over the free space value, there is 
presently no evidence for this. An upper bound for the 
global permeability of 1.8 is given by Behannon (1968). 
In the following discussion the values assumed for the 
global permeability and permittivity of the Moon are the 
free space values of these quantities. Possible departures 
of the permittivity from the free space value are insigni- 

j 

ficant in the subsequent application of this theory. 

First, consider the transverse magnetic (TM) 
or superscript e mode, characterized by the vanishing 
of the radial component of the magnetic field. The first 
step in obtaining the potential for this mode is to write 
out Maxwell's equations using (5) and (6) : 

VXE = k 2 H (7) 

VXH = -k^E (8) 
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where 


k 2 = icuti 


( 9 ) 


k = -k x k 2 


( 10 ) 


It can easily be seen that the radial component of (7) is 
satisfied exactly if 


E 


e _ 1 

<p r sin 0 


9U S 

3<p 


(ID 


e 1 auf. 

9 r d9 


( 12 ) 


where r, 9, <p are the spherical coordinate triplet and 
U e is a scalar potential. If (11) and (12) are next 
inserted into the 8 and c p components of (8) the 
resultant equations are 


a 

Sr 


K) = 


au 
w i ae 


( 13 ) 
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( 14 ) 


_i/ rH e] = zllL- in!. 

3r( r oj sinB 3<p 


By defining a potential Q e by 


U 


.A. JL 

kj. ^ 


(rO e ) 


(15) 


it can be seen that (13) and (14)' are satisfied exactly 
if 


„e _ 1 ao e 

w 8 sin0 3<p 


(16) 



an e 

30 


(17) 


Inserting (15) into (11) and (12) gives 


e = _L_ 5 2 (rfi e ) (18) 

E 0 k,r 3rd 0 


11 - 



( 19 ) 


G _X U , ^G v 

( p k^r s l n 9 Sr3<^5 r 

Finally, one must obtain an equation relating 
to fi e and a partial differential equation for Q e 
From the radial component of (8) together with (16) and 
(17) one obtains 


k E e 
1 r 


-v 2 n e 


_ 1 _ 

2 


_S_ 

Sr 


2 ML. 

Sr 


( 20 ) 


From either of the two transverse components of (8) 
together with (16) and (18) or (17) and (18) one obtains 



( 21 ) 


where 


k,E _ , Sk, .. . ,2 

lr ,2^e. 1 1 S , ^e x 1 S , _e % 

«=-F“- kn + 0 1F¥ (r “» -JTI 1101 (22) 

1 Sr 
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If separation of variables is assumed, then the only- 
solution for (21) is £ = 0 . With this result, 

equation (20) , and after some algebraic manipulation one 
obtains 



3k 


i a 


k-^r dr 3r 


+ -k 


(rQ e ) 


(23) 


and the determining partial differential equation for 
fi e : 

(V 2 + k 2 ) ^ ~(rCf) = 0 . . (24) 


It should be noted that the vanishing of the divergence 
of D (= eE) has not been invoked in order to obtain this 
set of solutions .7 This restriction was erroneously 
required in Wyatt's paper. In fact, the divergence 
equation is used to obtain the charge distribution which 
is not zero in the TM mode when the conductivity is a 
function of position. 

The complete solution of Maxwell's equation 
within the sphere requires a second set of solutions. 
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0 


the transverse electric (TE) or magnetic solutions, 
denoted by superscript m . In these solutions E™ = 
One obtains a set of potentials U™ and Q m in a 
manner completely analogous to the procedure followed 
above. In this case, however, the roles of (7) and (8) 
are interchanged. From the radial component of (8) one 
deduces a potential U™ such that 


m _ 1 SlT 

H e “ r S 0 


(25) 



l au* 

-r sin -0 d(p 


(26) 


Then using the transverse components of (7) with equations 
(25) and (26) one is led to the potential Q m such that 


* ‘ jq < 27 > 

f 9 - - is nra ' < 28 > 
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( 29 ) 


m _ 1 

E 0 sin 0 d<p 


r m = 1 S 2 

0 k.^r dr30 


(rn m ) 


( 30 ) 



k^r sin 0 dr b<p 


(rQ m ) 


( 31 ) 


Following the same technique as for the electrical 
potential one obtains for : 



+ k 


(rfi m ) 


( 32 ) 


and for the partial differential equation for 


0 


m 


(v 2 + k 2 )n ra = o 


( 33 ) 


The solution of equations (24) and (33) can now be 
obtained via separation of variables. The equations for 
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IQ 0 0 

determining ( E , K ) and (E , H ) from the potentials 
m 0 

of fi and Q are presented compactly in Appendix A. 

2.2.2 Spherical Harmonic Solutions in the Moon 

The solutions for equations (3) and (4) are 
easily obtained in terms of spherical harmonics by the 
method of separation of variables. The boundary condi- 
tions, used to fix the arbitrary constants of the solu- 
tions, are that the solutions are bounded at r = 0 and 
that the tangential components of the electric field are 
continuous at the lunar surface. For the transverse 
electric mode this is equivalent to the continuity of 
the normal component of the magnetic field (B_ = pH) . 

A discontinuity is permitted in the tangential magnetic 
field because of the confinement of the induced field by 
the solar wind pressure. This confinement is accomplished 
by a current sheet at the moon-plasma interface. 

rn 0 

The potentials Q and 0 are given by 

CD 

n m = pVH Q ^ sin <p P^A^G™ (r) (cos 0) (34) 


Q e = Hq ~ COS tp 3t^(r) pl t (co s 0) (35) 

1=1 
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where (r, 9, <p) are the spherical polar coordinates, a 

is the lunar radius, P^(cos 0) are the associated 
Legendre polynomials and 


_ i^(2Ul) 

“ t(Ul) ' 


i 




(36) 


m 0 

The functions G^ and are solutions of the differ- 

ential equations 


d 2 G^ 

dr 2 


k" - 


1 ( 1 + 1 ) 

2 


K (r) 


= 0 


(37) 


d M 


dr 


4^!^l + / k 2 .iiM, 

k 2 dr dr I r 2 j 


G e t (r) 



(38) 


0 

In Sonett et. al. (1971a) the appearance of A^ and A^ 
were suppressed, these constants being incorporated into 
the definitions for G^ and . . The exact form of 

the boundary condition at r = a is obtained by using 

the spherical harmonic expansion for the incident 

/ 

and internal fields, given in Appendix B. The boundary 
conditions for the fields are given by 
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(39) 


A™ G™(a) = j t (2va/-\) 


= ^T" A ( a i, (2iraA)) (40) 

4tt a da ^ 

r=a 

where j^(x) is the spherical Bessel function of order 
l . 

m s 

2.2.3 Behavior of the Solutions G^ and Near 

r =.0 

In the center of the Moon, near r = 0 , it 

dk 

is reasonable to assume that ^ becomes vanishingly 

m g 

small. Then G^ and G^ in the vicinity of r = 0 
are of the form 


dGT 
e l 

dr 


T • _ m / g / » 

Lim G. (r) 

r-»0 


r 3<, (k o r) 


(41) 


where k Q is the propagation constant at r = 0 
Therefore, at r = 0 


G® ,e (o) = 0 


(42) 
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(43) 

For numerical calculations the boundary conditions 
expressed by (42) and (43) are useless. However, an 
auxilliary function can be defined which makes the solu- 
tion of (37) or (38) more amenable to digital computer 
techniques. Define the function N^' (r) by 



G™' e (r) = r U yj' e (r) 


(44) 


Then from (37) and (38) the differential equations for 
XT m,e 

are 


-ai ^ + Hitii p + k y? ( r> . 0 

^ r 2 l r dr t 


(45) 


j2 e 

Zl J 2(1+1) 

, 2 + I r 
dr 


d . 2. 

5? (lnk > 


dr 


1 


k2 - r ^ (lnk2) (’ N ! (r > =0 


(46) 
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and by using (41) the boundary conditions at r = 
become 


<' e ( o) 


2 l llk l 
o_ 

(2i+i ) : 



m, e 
l 


r=0 


0 


m 0 

The constant in (47) can also be absorbed in A^' 
give as an alternate condition to (47) , 


N 


m, e 
l 


(o) 


1 


The boundary conditions at r = a then become 


A™ = j ^(2iraA)/ 


t+1 m . 
a N^(a) 


( 47 ) 


(48) 


to 


(47a) 


(49) 
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A; = 


k 2 X 2 /4TT 2 a]^-[aj^ ( 27 ra/X) 


dr 


|r t+1 N J(r)| 


(50) 


r=a 


For those cases where the differential equa- 
tion (either (45) or (46) ) can be integrated over the 
entire range 0 ^ r ^ a the use of the auxiliary func- 
tions provides a means for reducing the problem to the 
solution of an initial value computation. Care must be 
taken, however, because of the nature of the differential 
equation. These equations are basically diffusion 
equations , the solutions -becoming -exponentially damped 
as r decreases from a toward zero. For higher 
frequencies and, hence, smaller local skin depths 
6 , the attenuation may become so great that 

significance will be lost in the integration scheme. 

For these cases an alternative procedure would be to 
integrate the original differential equation ((37) or 
(38)) from some depth b where the local skin depths 
6 is small (6 ^ (a-b) ) , and assume that the fields 

vanish for r < b . The equations can then be solved 
as a boundary value problem with N^' (b) = 0 and 

N^' m (a) determined through equations (49) and (50). 
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3.0 


The Near Surface Lunar Thermal Gradient 


3 . 1 Introduction 

The initial results from the Apollo 12 lunar . 
surface magnetometer experiment show that on the sunward 
hemisphere the moon displays a magnetic field response to 
the interplanetary electromagnetic field which is larger 
by one-half order of magnitude than an equivalent vacuum 
bounded case. This strong interaction provides the basic 
justification for the current layer at the moon-plasma 
interface, which is used in the theoretical development 
in Section 2 . Without this layer there can be no signal 
increase over the vacuum case. The field increase or 
amplification is dependent on the combined exclusion of 
the induced field by a highly conducting lunar core and 
the confinement at the sunward surface of this, same field 
by the solar wind pressure. 

The interior field exclusion is related to the 
damping of the TE wave as it propagates into the moon. 

For a constant conductivity moon the damping occurs at 
progressively shallower depths as the frequency increases 
This phenomenon is accentuated when the electrical con- 
ductivity increases with depth. Even for an inverted 
conductivity profile one can easily be convinced that 
the magnetic field compression at the lunar surface in- 
creases with increasing frequency. 
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The first LSM data points for the magnetic field 
magnification were obtained for frequencies at the high 
end of the usable band for the Explorer 35 Ames magnetometer 
and LSM combination , .02 to. 04 Hertz. .• These frequencies are 
high enough to lead to the expectation that the electromag- 
netic damping in the lunar crust will occur at some depth 
d where the lunar temperature profile can be described by 
a linear thermal gradient approximation; i.e., 

T(r) = T g + T r (a-r) , b < r < a (51) 

where b is the smallest radius for the approximation and 
where T is the surface temperature in degrees Kelvin, 
a is the lunar radius (1740 kilometers), T r is the thermal 
gradient in °K/km and r is the radial distance in kilo- 
meters. The estimates of the thermal gradient are depen- 
dent on the model conductivity assumed although there are 
obvious constraints provided by the range of validity of 
the linear gradient approximation and necessary high 
internal temperatures for the low conductivity models. 

Since only the electrical conductivity a 
occurs in the differential equations for the field 
(equations 37 and 38) rather than T(r) the model depen- 
dence is obvious. As a final point, one should note that 
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the magnetic field exclusion occurs at depths where the 
local skin depth, 6 = ~\j2/udfro , is of the order Of one 
kilometer and the quantity (a/ cue) is much greater than 
unity. Hence, the equations (37) and (38) are basically 
diffusion equations and the value of e is of no conse- 
quence to the resultant field distribution. By combining 
equation (51) with various models for rock electrical 
conductivity and then comparing the resultant calculated 
magnetic field with the experimental data, a series of 
estimates can be obtained for the lunar thermal gradient 
near the surface. 

3 . 2 Theory 

The refined data from the LSM experiment is 
generally present as a ratio of the field measured at the 
LSM to the inducing field measured at Explorer 35 . This 

■ V 

data is further obtained only after Fourier decomposition 
of both original signals. Thus, the experimental data 
provides a transfer function of the moon to electromag- 
netic induction. In order to make a comparison between 
the theoretical and experimental data, a theoretical 
transfer function must be defined. This is easily accom- 
plished for the TE mode by going to Appendices A and B 
of Section 2 and equation (34) . First, remember that the 
boundary condition for the TE mode requires that the 
radial components of the magnetic field be continuous. 
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Then the transfer function must relate the tangential 
components of the field. A normalized harmonic transfer 
function is defined for the TE mode by' 



(r) 


d 

dr 


Id,. 

a 


m ! \ 
G t (f> 

^(2-Tra/A) ) 


( 52 ) 


This provides the ratio, for each harmonic, of the TE 
mode magnetic field anywhere in the moon to the inducing 
magnetic field in that same harmonic at the lunar surface. 

This definition for a transfer function which is indepen- 

/ 

dent of 9 and ,<p cannot be made for the TM mode mag- 
netic fields. Instead, the induced TM magnetic fields 
can only be compared to the inducing TE magnetic fields. 
The difficulty that arises can be seen by examining the 
limiting case as the frequency goes to zero. For this 
limit the magnetic field components of the driving field 
for the TM mode vanish while the induced TM components can 
have any value (theoretically) depending on the electrical 
conductivity of the moon model. This would lead to 
extremely large values of a transfer function even though 
both the driving field and induced field are well behaved. 
As an alternative, the TM transfer function has been 
defined as 
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( 53 ) 




2-rria^ 

Tvr 


G t^) 

^(aj t (2ira/X>) 


. e 

For the models used in this section, the values are 

too small to be of importance. 


3 . 3 Conductivity Models 


The electrical conductivity of rocky material 
is conventionally represented by a function of the form 


-I *0i 

i=l 


exp(- ^) 


(54) 


where is the mobility (conductivity for T = 00 ) , 

is an activation energy, the energy required to raise 
an electron from a valance band into a conduction band, 
k is Boltzmann's constant, and T the absolute temper- 
ature. The dependence of conductivity upon temperature 
and activation energy as expressed by equation (54) 
follows the classical theory of semiconductors, but 
information so far developed in geology lags the former 
significantly. Nevertheless, from an experimental stand- 
point, equation (54) does represent conduction sufficiently 
well so that application to the present problem can be 
made. The total conduction is made up of an impurity 
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range, one in which intrinsic electronic conduction 
dominates and lastly an upper temperature range in which 
ionic conduction takes over. The latter is likely not 
significant in the application of this paper and is dis- 
regarded. The use of a conductivity function at atmos- 
pheric pressure extrapolated to depth can have the effect 
of reducing the apparent conductivity because' of the 
presence of interstitial cracks. We have no information 
on this for lunar samples. Any error would be in the 
direction of increasing the calculated thermal gradient 
so that this must be regarded as a limitation, as is the 
paucity of samples and the generalization to bedrock. 

We consider here two conductivity functions 
thought to represent extremes for the lunar material. 

For an Apollo 11 lunar crystalline rock (sample 10024-22) , 
Nagata et. al., (1970) found 

Gj/k = 5800 °K , e 2 /k = 14,500°K , 

g 

<Jq^ = 7.9 mhos/m and cr^ = 3.7x10 mhos/m . It has 
been found that an Apollo 12 rock (sample 12053-47) has 
a. more complicated dependence of conductivity upon tem- 
perature but displays a generally similar behavior to 
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that of the Apollo 11 sample with somewhat reduced con- 
ductivity over the range measured. 

There is good reason to suspect that these 

samples are not representative of the basement rock. 

For example, they could not constitute the basic rock of 

the Moon for the basalt transforms to eclogite at 12 kb 

pressure, and the density of this polymorph is too high 

to match the pverall lunar density. Other more subtle 

geochemical arguments having to do, for example, with the 

incompatible element concentration in the basalts, suggest 

that the interior is better represented by a periodtite 

containing a substantial fraction of olivine. For this 

reason we have chosen as the other end point calculation 

an olivine conductivity function (England et. al . , 1968) 

for which e-j/k = 10,670°K , e 2 /k = 31,300 °k , = 

7 

55 mhos/m and ct Q 2 = 4x10 mhos/m . The behavior of 
the Moon rock and olivine as functions of temperature 
are illustrated in Figure 1. 

3 .4 Numerical Calculations 

For sufficiently small skin depths it is 
unnecessary to integrate the equations of the previous 
section to the center of the Moon. It is sufficient to 
replace the deep interior of the Moon by a perfectly con- 
ducting sphere and integrate the equations from the lunar 
surface to the surface of the perfectly conducting core. 


The exact location of this lower boundary is unimportant 
as long as the local skin depth is small . For convenience 
in the numerical, computations to follow, the lower boundary 
has been set at a depth of 250 km. The transfer function 
was obtained for several values of the frequency, f , and 
a solar wind speed of 400 km/sec . 

The equations derived in Section 2 have been 

numerically integrated for the lunar thermal gradient 

Moon models described above. For all models T was 

s 

r 

taken as -30 °C or 243 °K. The behavior of the real and 
imaginary parts of the transfer function for the -1 = 1 
TE mode (equation (52)) is shown in Figure 2 for the 
Nagata-Moon rock conductivity, for a frequency of .04 Hz 
and thermal gradients of 2, 3 and 4°K per kilometer. 

Except. „f or the lowest thermal gradient, the magnetic 
fields and, hence, the transfer function are well damped 
below r = 1500 kilometers . 

Extensive calculations for the TE mode, 1 = 1 , 
transfer function have been made for the two conductivity 
functions. The data from these calculations is shown in 
Tables 1 and 2 where the magnitude and phase of the 
transfer functions are given for frequencies, of .02, .03 

and .04 hertz and for thermal gradients from 1 °K to 7°K 
per kilometer for the England-olivine conductivity. Data 
was calculated only up to 4°K/km for the Nagata-Moon rock 
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i 

1 

I ' 

| T ( °K/km) 
i r 

Freq (Hertz) 

1 

.02 

i 

.03 

.04 

l 

1 

2 . 4/- 9.5° 

2 . 6/- 9 . 8 ° 

2 . 9/-10 ° 

2 

3 . 9/ - 9 0 

4 . 4/ -9.2° 

4 . 7/ -9.7° 

3 

5 . 3/- 8 .3°. 

6 . 1/-8 . 6 ° 

6 . 4/- 9.0° 

■ 

4 

6 . 7/ -8 . 0 ° 

7 . 5/-8 . 3 ° 

r 

8 . 1/-8 . 4 ° 


TABLE 1 

TE Mode t = 1 transfer function 
calculated using the Nagata-lunar 
basalt conductivity 
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Freq (Hertz) 


• 

T r (°K/km) 

.02 

.03 • 

.04 

1 

1 . 5/- 6 . 1 ° 

1 . 6/- 6 . 2 ° 

1 . 7/-6 . 2 ° 

2 

2 . 4/ -7.8° 

2 . 6/- 7.7° 

2 . 8/-7 . 7 ° 

3 

3 . 3/-8 . 3 ° 

3.6/-8 .2° 

3 . 9/ -8 . 2 ° 

4 

4 . 1/-8 . 5 ° 

4 . 5/-8 . 5 ° 

4 . 8/ -8.5° 

5 

4 . 9/-8 . 8 ° 

5 . 3/-8 .7° 

5 . 7/-8 . 6 ° 

6 

5 . 7/-8 . 9 ° 

6 .-2/-8 . 8 0 

6. 6/- 8 .7° 

7 

6.4/-9.0 0 

7.0/-8 .9° 

7 . 4/-8 . 9 ° 


TABLE 2 


TE Mode,. t = 1 , transfer function 
calculated using the England-olivine 
conductivity 











because of its much higher conductivity and resultant 
large transfer functions . The amplitude data of 

Tables 1 and 2 are displayed graphically in Figure 3. 

From this figure it is clear that for a given value of 
T™ at a given frequency the olivine conductivity requires 
a higher thermal gradient and, hence, a higher internal 
temperature. Calculations were also carried out for 
higher order harmonics for the TE mode. These results 
listed in Table 3 were obtained for f = . 04 Hertz, a 

5 

solar wind velocity V of 4x10 m/ s and a lunar radius 
0 

of 1.74x10 meters. For these values the argument of 
the Bessel functions in the driving term (2irfa/v = 2-jra/A) 
is approximately unity. It can be seen that the ampli- 
tude of the transfer functions and, hence, the amplifi- 
cation of the field decreases with increasing harmonics. 
In addition, for this set of parameters, the higher order 
harmonics in the forcing functions are of decreasing 
importance in the harmonic expansion of the forcing field 
at the lunar surface. These forcing functions are in the 
ratio 1 . 0 :0 . 18 : 0 . 025 : 2 . 8xl0 - ^ : 2 . 5xl0 - ^ . Only the second 
harmonic produces any appreciable correction in the field 
The magnitude of this correction will decrease with 
decreasing frequency since for 

r • • , , . (2x) n n! 

3 n (x) - (2n+1) . 

X -* o 
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England-Olivine 



TABLE 3 

Normalized Surface Tangential Transfer functions for 

f = .04 Hertz, v = 4x10^ m/s,, a = 1.74x10^ meters. 

In the harmonic expansion the forcing functions for 

-3 

the harmonic are in the ratio 1:0.18:0.025:2.8x10 : 2. 5x10 
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3.5 


Discussion 


The response functions shown in Figure 2 have 
the general form one would expect for magnetic induction 
in a layered spherical conductor for which the inner con- 
ductivity is much greater than the conductivity near the 
surface. The phase change with increasing depth is a 
characteristic of the change of the governing differential 
equation from one which is essentially a Laplace equation to 
a diffusion equation. The strong damping of the induced 
field starts in the region of the phase change. It is 
clear that the field penetration decreases as the thermal 
gradient increases. The property that the induction is 
relatively constant "in the outermost regions corresponds 
to the fact that in this region both the real current 
J = oE and the displacement current edE/dt are negli- 
gible. The results of the calculation are qualitatively 
similar to results which would be obtained for a two 
layer calculation. For instance, if the Moon had infinite 
conductivity for r < b and zero conductivity for 
b < r < a , the transfer function T m (r) for those fre- 
quencies where |27ra/A|« 1 is given by 


T m (r) = L + l/2Lb/rjl b < r < a 
1 - <b/a) J 


(55) 
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From this equation it is obvious that if b is close 
to a , T m (r) is very nearly constant and equal to 
1 . 5/ (1 - (b/a) 3 ) . , 

The magnitude of the transfer functions given 
in Tables 1 and 2 and shown in Figures 2 and 3 are indica- 
tive of the compression of the induced magnetic field 
between the solar wind at the lunar surface and the highly 
conducting lunar interior. For a given conductivity 
function and frequency the field amplification increasing 
with increasing thermal gradient. Similarly, for fixed 
thermal gradient and conductivity the magnification 
increases with frequency. 

The penetration of the magnetic field into the 
Moon increases as the frequency decreases for a fixed 
thermal gradient and conductivity function or as the 
thermal gradient decreases for a fixed frequency and con- 
ductivity function. The decrease in the magnetic field 
amplification corresponds to the increase in penetration. 
Since the true thermal gradients must be expected to 
decrease with increasing depth, the lower .frequency waves 
will penetrate even further into the Moon than is indi- 
cated by the constant gradient analysis. This analysis 
is only meaningful for the highest frequencies where the 
major damping occurs within the first 200 to 300 kilo- 
meters. The depth of penetration has been obtained for 
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TABLE 4 

Depth of penetration of the induced field 
using the first zero of the real part of 
the TE , t = 1 transfer function 
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It can be seen from this table and from Tables 1 and 2 
that similar values of the damping depth give approxi- 
mately the same value for the transfer function. 

3.6 Comparison with Experimental Data 

The presently accepted experimental values for 
T™ (a) are shown in Figure 5 with the error bars (Sonett 
et al. 1971a and 1971b). The experimental values for .02, 
^03 and .04 Hertz are shown in the curves in Figure 3 with 
appropriate error bars. An immediate difficulty arises 
because the value of T™ for .02 Hertz is greater than for 
the other two frequencies. There is, as yet, no theoretical 
explanation for this. However, using the data for .03 and 
.04 Hertz a lower limit can be obtained for the thermal 
gradient of 1.4°K/km. This is obtained for the Negata- 
basait. If one uses a less conducting material such as 
the olivine, a higher value of the thermal gradient is 
obtained, 2.7°K/km. 

It has been suggested that some TM mode 
admixture may be occurring. To investigate this possibility 
the TM mode £ = 1 transfer function was calculated for a 
frequency of .03 Hertz and thermal gradient of 4°K/km 

-4 

using the Moon basalt conductivity. A value of 2.5x10 /0° 

was obtained for the transfer function showing that the 
low conductivity of the crust was cutting off the 
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excitation of this mode within the Moon. To the conduc- 


tivity equation given by (54) with the basalt values a 
constant conductivity was added; i.e., 


ct(T) 


Basalt 


(T) .+ o 


(56) 


The results of these calculations are shown in Table 5. 


TABLE 5 

TM transfer functions for f = .03 Hertz ,! 

j 

T r = 4°K/km and Nagata Basalt 


• 

. 

m e 

a o 

T 

0 

2 . 5xl0~ 4 /0 ° 

CO 

1 

o 

p — i 

7 . 7xl0~ 3 / 0 ° 

[•'- 

1 

o 

•H 

.072/-. 17° 

io " 6 

. 56/- 1 . 25 ° 


It can be seen that the TM mode will not produce any 
measurable effect unless the surface conductivity is at 
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least of the order of 10 mhos/meter. This corresponds 
to a temperature of 365 °K for the Moon basalt and 590 °K 
for the olivine, much higher than the mean surface lunar 
temperature of 243 °K. 



4 . 0 The Conductivity of the Moon 

4.1 Introduction 

In the previous section an attempt was made to 

bracket the near surface lunar thermal gradient. In this 
section the total spectral data obtained from the LSM 
experiment is used to provide an estimate for the conduc- 
tivity of the moon from approximately r = 400 km to the 
surface. The general properties of both the LSM and Ames 
magnetometer on Explorer 35 have been given (Dyal and 
Parkin, 1971, Milhalov et. al, 1968) . The theory has been 
described in Sections 2.0 and 3.0. This work is also 
included in Sonett et. al. (1971b). 

4 . 2 Experimental Data 

The spectral determination and coordinate trans- 
formations for the experimental data were carried out at 
Ames Research Center (Sonett et. al. 1971b) . The forcing 
field defined by Explorer 35 measurements is transformed 
into the local LSM coordinate system (x is along the normal 
outwards from the surface, y is easterly and z is northerly 
at the site of ALSEP and the LSM) . Power spectral densities 
are obtained for the magnetic field time series at both 
Explorer 35 and LSM for all three components. 

A representative set of spectra at both Explorer 

and LSM for the y and z magnetic field components (tangent 

to the surface) is shown in Figure 4. The forcing spectra 

-2 . 

are seen to display the expected f dependence characteristic 
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•of interplanetary fluctuations, while the LSM shows an 
f dependence. Thus, significant frequency dependent 

l 

power amplification is apparent in the record. A combi- 
nation of all the data points A(f) =ji (A^(f) + A^(f))j ^ 
for all frequencies up to 0.040 Hz has been made. See 
Figure 5 . 

4 . 3 Lunar Electrical Conductivity Profile 

The theory of lunar induction, discussed earlier, 

has been used to derive a profile of electrical conduc- 

( 1 ) 

tivity from the empirical transfer function. The theo- 

retical amplification for the lowest TE mode, computed by 
numerically integrating equation (45) with a radially de- 
pendent conductivity, is matched to the empirical trans- 
fer function using a Newton-Raphson iterative scheme. 

This scheme readjusts the entire conductivity profile to 
yield amplifications which best fit the data in a least- 
squares sense over the entire frequency range. This has 
been carried out using frequency values of 0.83, 1.75, 5, 

12, 17, 22, 25, and 35 millihz. The conductivity profile 
is characterized by its values at the radial locations 
r = 800, 1200, 1400, 1450, 1490, 1510, 1550 and 1740 km. 

For r <800 km the conductivity d is set to the value at 
r = 800 km; elsewhere a linear interpolation of log o is used 

(1) This theory, and the computer program for deriving 
• the conductivity profile from the empirical transfer 
function were developed by American Nucleonics Corp. 
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The computer calculation is started with a 
continuous conductivity profile defined by 8 parameters, 
as described above. Numerical integrations are carried 
out to obtain values of amplification at the eight fre- 
quencies. A comparison of these amplifications with the 
data provides the Newton-Raphson scheme with the input 
required to adjust the 8 conductivity parameters to yield 
a profile whose corresponding theoretical amplification 
curve is a better least-squares fit to the empirical 
transfer function. The iteration is continued until an 
adequate fit is obtained. The result of such a fit to 

— f 2. 2 2 1 

the empirical amplification curve A = < — (A^ + A^)) 
is shown .in Figure 5 . The dif ferences between the theo- 
retical amplification and the empirical A are very 
likely attributable to the various complications of the 
excitation process not accounted for by the theory, as 
well as considerations of computer running time which 
limits the number of frequencies and conductivity para- 
meters (spatial resolution) employed. 

A conductivity profile is associated with the 

theoretical amplification curve of Figure 5 and is shown 

in Figure 6. The conductivity profiles derived from the 

iterative least squares inversion for A , A and their 

y z 

on standard deviation limits (see Sonett et. al., 1971, a 
orb). The prominent spike is an invariant characteristic 
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of the inversion for each set of data. It is centered 

at about r = 1500 km where the conductivity is nearly 
-2 ' 

10 mhos/m. The inner minimum lies at about r = 1400 km 
and the conductivity appears to rise at greater depth. 

Gross bounds on the conductivity profile are seen in the 
insert; these are determined from the one standard devia- 
tion limits of the various A's, but do not themselves 
represent one standard deviation limits on the conductivity 
profile. 

The computer calculations for Figures 5 and 6 

-4 

started with a constant conductivity of 10 mhos/m 
However, a number of computations have been carried out 
using different values for uniform starting conductivities 
(e.g., 10 and 10 mhos/m) and different radial 
locations. In every case tested, initial convergence 
was rapid and the final conductivity profile invariably 
displayed the prominent spike near r = 1500 km . For 
the profiles reported here, several values of r were 
chosen in the neighborhood of r = 1500 km to better 
define the conductivity spike. 

Whereas the large spike in conductivity is a 
persistent feature of the inversions, the character of 
the conductivity profile at greater depth, where the 
conductivity appears to rise, is not so certain and our 
results for the conductivity at these depths must remain 
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tentative. The surface amplification is a rather insensi- 
tive function of core conductivity because of the relatively 
small core volume and the distance to the surface. 

4 . 4 Comparison with Other Conductivity Profiles 

A physical understanding of the appearance of 
the large conductivity spike can be obtained as follows. 

In a two layer model with an infinitely conducting core 
and a non-conducting shell the amplification of the tan- 
gential magnetic field components' equation (55), can be 

written in the form 

_ 1 3 ( Core volume) 

2 (Shell volume) 

This simple result requires that the high frequency skin 

depth barrier be near r = 1500 km for an amplification 
« 4 . At frequencies above 0.02 Hz essentially no wave 

penetration takes place through this barrier. Thus, the 
volume available in the core for field lines is insigni- 
ficant compared to that in the non-conducting shell. 

This means that the amplification becomes independent of 
frequency; i.e., dA/df ^ 0 at high frequency. Con- 
current with the high frequency limitation is the require- 
ment that A drop to near unity at the lower frequencies. 

If it were assumed that the electrical conductivity were 
mono tonic, so that the interior conductivity were uniformly 
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high, then the low frequency amplification would be in 
excess of the observed values. 

The inadequacies of monotonic conductivity 

profiles are shown in Figure 7 . The experimental r.m.s. 

amplification data A (circled points with associated 

error bars) can be compared with the amplification curve 

for the A conductivity profile of Figure 6 (labelled 

"best fit") and amplification curves of several monotonic 

conductivity profiles. The "2. layer" amplification curve 

was calculated for a model with a core of radius of 1560 km 

-4 % 

and a constant conductivity of 7.6x10 mhos/m , and a 

shell of zero conductivity. This value of the core 

conductivity provides a best fit to the experimental A 

for the given core radius. Other best fit two layer 

models with different core radii have been investigated; 

the one used in Figure 7 yields an amplification curve 

which matches the data most closely. The amplification 

curve labelled DYAL and PARKIN (1971) was computed from 

a 3 layer conductivity model proposed by these authors. 

The parameters of their model are a = 10 mhos/m for 

-4 

0 < r < 1044 km , a = 1.7x10 mhos/m for 1044 km < 
r < 1653 km and a = 0 in the outer shell. In other 
3 layer conductivity models consistent with the data 
analysis of DYAL and PARKIN (1971) , the outer boundary 
of the region of intermediate conductivity can be located 
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anywhere between r = 1653 and 1740 km . We have 
computed the amplification curves for a number of these, 
additional 3 layer models; none provide a better fit to 
the observational curve. Also shown in Figure 7 are 
amplification curves from conductivity profiles proposed 
by SILL (1971) (model 3 in that paper) and NESS (1969) . 
Other conductivity models investigated by SILL (1971) fit 
the experimental data no better than his model 3. The 
conductivity model of NESS (1969) consists of a core of 
radius 1426 km with a = 8x10 5 mhos/m , and a non- 
conducting shell. Amplification curves from conductivity 
models of WARD (1969) were also considered but these were 
rejected upon comparison with the observational data. 

The "2 layer" amplification curve shown in 
Figure 7 is, compared to the experimental data, high at 
f < 0.0035 Hz and f > 0.03 Hz , and low in the range 
0.01 Hz < f < 0.02 Hz . The data show a flat response 
at frequencies above 0.02 Hz, whereas the slope of the 
"2 layer" curve is high at these frequencies. Further- 
more, as the frequency increases the curvature of the 
data changes from positive to negative at about 0.003 Hz. 
The "2 layer" amplification curve is everywhere concave 
up. This 11 s" shaped character of the experimental data 
is faithfully reproduced only by our "best fit" conduc - 
tivity profile . Further work is in progress to improve 
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the high frequency slope of our "best fit" model. This 
will probably lead both to a higher value of the maximum 
conductivity and a steeper slope of the conductivity 
profile at the outer edge, effects which will tend to 
sharpen the conductivity spike in order to allow pene- 
tration of low frequency magnetic field fluctuations. 

4 . 5 Compositional and Thermal Model 

The electrical conductivity profile cannot be 
explained by a uniform material and a plausible thermal 
profile. In the region from the surface to r = 1500 km 
where the conductivity attains its maximum value, the 
rise of conductivity with depth is a reasonable consequence 
of the accompanying increase of temperature in a material 
of uniform composition. Below R = 1400 km the apparent 
rise in conductivity is again explainable by an increase 
in temperature. On the other hand, the precipitous 
decrease of electrical conductivity by 2 to 3 orders 
between r = 1500 and 1400 km cannot be explained as due 
to temperature. Either a compositional change, phase 
change, or a combination of the two is required. Thus, 
a reasonable model for stratification of the Moon, limited 
by the present poor spatial resolution of the analysis 
indicates a core out to r « 1400 km overlain by a mantle 
of higher conductivity material, plus possibly a transition 
layer at r = 1400 - 1500 km . 
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In order to infer a temperature profile from 
the conductivity profile it is necessary to use conduc- 
tivity-temperature functions of known rock materials. 

The Apollo basalts (Schwerer et. al., 1970; Nagata et. al 
1970) are representative of the most conducting rock 
matter known. On the other hand, olivine (dunite) or 
olivine-peridotite (England et. al./ 1968) is represen- 
tative of poorly conducting geological material. The 
peak conductivity value found in the lunar mantle corre- 
sponds to a temperature of about 450 °C for lunar basalt 
or 950 °C for olivine. Since the temperature just under 
the lunar surface is -30 °C, the corresponding thermal 
gradients in the outer mantle are .about .2 °K/km and 4°K/km 
These estimates correspond to those obtained in Section 
3.0 for .02 Hertz. 



5.0 


Early Planetary Heating by Electromagnetic 


Induction 

5 . 1 Introduction . 

Sonett et. al. (1968) and Sonett et. al. (1970) 
have given convincing argument for the possibility of a 
short, intense heating episode caused, by electromagnetic 
induction, this episode occurring during the initial 
formative years of the planetary body. In these two 
papers the electrical heat source was the unipolar induc- 
tion field which is the d.c. limit of the TM mode discussed 
in Section 2.0. The model parameters have been thoroughly 
presented in the two references and in Schwartz (1969) . 

An obvious deficiency of the theory, as presented to date, 
is the reliance on only the d.c. field component to the 
exclusion of the rest of the frequency spectrum. To 
remedy this situation, the original planetary heating 
computer program has been rewritten to provide a method 
for evaluating the separate and combined contributions 
to the planetary heating from the TM and TE modes as a 
function of frequency. The program is still in the check- 
out stage as far as practical results are concerned. The 
individual subroutines have all been checked out but the 
actual modes of operation have not been worked out. Some 
brief runs have been made to ensure that all the routines 
work properly together. 
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5.2 


The Driving Field 


As described in Sonett et. al. (1970) it is 
assumed that the sun passed through a T Tauri stage 4.5 
billion years ago on the way to becoming a main sequence 
star. During T Tauri the sun maintained both a high spin 
rate, near the centrifugal limit, and a high surface 
magnetic field. The driving field for the electromagnetic 
induction is assumed to be the quasi steady state rotating 
sector structure of the magnetic field and the frozen-in 
V x B electrical field which is convected past the planet 
by the solar wind. A two level sector structure for the 
magnetic field is assumed so that the magnetic field H 
and also the electric field E = V x pH have magnitude 
H q and E q , respectively, but change sign (direction) 
every half period of the solar rotation. If the solar 
period is 2 t(t = Tr/a> 0 where a> Q is the local time solar 
spin rate) then the time dependence of driving fields is 
given by 



2nT ^ t £ (2n+l) T 
(2n+l) T ^ t £ 2 (n+1 ) T 


(57) 


- 50 - . 



The average power deposited by joule heating can be 
obtained using Fourier Series techniques with the result 
that for any given mode (i.e., TE or TM) 


H(r) 


+ to 



n=-a= 


C V (r) 
n n 


2 


(58) 


3 

where H(r) is the power density in watts/m , cr(r) is 
the electric conductivity in mhos/meter, V n (r) is the 
normalized electric field in volts/meter and C n is the 
Fourier coefficient for the driving term 


C = 
n 


/ E ) 



< 

\ H 


\ Ol 



0 n even 
-2i 


mr 


n odd 


(59) 


The square of the total input electric field is 


09 

= 


(2m+l) 2 = E 2 
o 


m=o 


(60) 


- 51 - 


Thus, the total driving field looks, in magnitude, iden- 
tical to the driving field for a constant field of E 

o 

In fact, the total input for the TM portion of the driving 
field is equivalent to the input for the unipolar induc- 
tion where one assumes C =1 and C =0 for n = 0 

o n 

The highest solar spin rate, at t = 0 years, 

-4 

is (go ) = 4.2x10 radians/sec or (f ) = 

o max c max 

6.67x10 Hertz . If one sums the magnitude from the 
first 10 terms of the TM mode, the omitted terms consti- 
tute 2 % of the total used. 

5 • 3 General Program Concept 

The main program, called HEAT03 ■, solves the 
basic problem of the thermal history of a planetary body. 
Input data is read into the program in the main body and 
through a subroutine CONIN. Provision is made for 
radio nuclide heat sources in the subroutine set RADIO 
and AHEAT. The electromagnetic heat sources are obtained 
through the subroutine ERAD and a subsidiary subroutine 
FSET. The subroutine ERAD sets up the T Tauri flow 
parameters, the driving electromagnetic field, solar spin 
rate and the electrical conductivity. The calculations 
of the electric fields and the joule heating are carried 
out through FSET. In the present version of the code, 
the calculation can be made for: 
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(1) the d.c. unipolar generator 

(2) the TE mode for a series of n frequencies 

(n s 1) where f = (2n-l) f where f 

n o o 

is the solar spin rate 

(3) the TM mode for the same series of n 
frequencies as (2) . 

Combinations of (1) and (2) or (2)' and (3) can also be 
obtained. For the present n is set equal to 1 in 
FSET because more investigation is needed to determine an 
optimum method for running the program. The program is 
designed for no melting and therefore no radioactive 
concentration. The radioactive sources are assumed to be 
uniformly distributed within the body. The electrical 
heating term has two parts, one from the tangential 
electric field and one from the radial electric field 
(TM mode only) . The program can be set through input 
data to use either or both parts. 

5.4 Input Deck 

The following is the make up of a typical run. 
Each card set is explained. 

Card 1 This is a descriptive header card using columns 
2 through 80. It must be used but can be blank. 

Card 2 Name List Ml. There are three inputs on this 
name list. 
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ATIM is the start time of the calculation 
in years . 

TTESTZ is a skip time in years. For TTESTZ 
= 6 > 0 the electric field joule 
heating will not be calculated any 
more often than 6 years, except 
when a printout of the thermal profile 
and the heat sources is called for. 

EFCT determines which components of the 
electric field will be used for the 
joule heating. 

1. E , oE 2 

r r 

2 . E , • i , oE 2 

Tangential T 

J 

3. a(E 2 + E 2 ) 

r T 

Card Set 3 Name List M2. This name list has 15 inputs. 

GTC This is the time constant, y, for 

the solar magnetic field in units of 
-1 9 

eons (1 eon = 10 years) . The 
maximum permissible value is .457. 

For any larger value the sun becomes 
centrifugally unstable at time zero 
according to the despin theory used 
for the program. The magnetic field 
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has the variation 


RBAR 


RG 


DP 


ENZ 


B ( t) = B 4 5 exp y (4 . 5-t) 

This is a parameter which might just 

as well be an internal constant, equal 

2 8 
to r /r where r = 6.99765 x 10 
s se s 

meters solar radius. r = 1.495 

se 

x 10^ meters - earths distance from 
sun. RBAR = 3.2754 x 10^ meters. 

This is the distance between the sun 

» 

and the planetary body under investi- 
gation, measured in astronomical units 
(a.u.) . For the earth, RG = 1. 

This is a constant term which is 
added to the electrical conductivity 
in ERAD. DP can have any value (in 
mhos/meter) but for the most recent cal- 
culations, DP = 0. 

This is the particle flux density of 
the solar wind at 1 a.u. and at 
t = 4.5 eons. The present value is 

7 

taken as 10 protons per cubic meter. 


-5 5- 



ENTTZ 


This is the particle flux density for 


GSW 


GTT 


VTZ 


VELZ 


GVZ 


TS 


TTS 


GTTS 


a T Tauri flow at 1 a.u., t = 0 , 

and normalized to a solar mass loss 

-30 

Am = M q = 2 x 10 kilograms . For 

the assumed velocity and time decay, 

14 

ENTTZ = 8.6 x .10 protons per cubic 
me ter. 

Time constant for the solar wind flux 
density, taken as 1.0214 per eon. 

Time constant for the T Tauri flux 
density, taken as 1200 per eon. 

T Tauri flow velocity at t = 0 : 

5 

2 x 10 meter/sec. 

Solar wind flow velocity at t = 

4.5 eons; 4 x 10 5 meters/sec. 

Solar wind velocity time constant; 
85.5/eon. 

Average space temperature at 
T = 4.5 eons - normal, 273° Kelvin. 

Increment in space temperature caused 
by T Tauri obscuration; 0 ^ TTS ^ 500 

Time constant for decay of obscuration 
1.2 x 10 per year. 
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Card Set 4 Name List M3 has six entries. 

DELA Spatial difference step for inner 

portion of the planetary body in the 
solution of the thermal diffusion 
equation/ normally 2 x 10 meters . 

DELB Spatial difference step for the outer 

layer of the planetary body in the 

solution of the thermal diffusion 

equation. For early heating with T 

5 

Tauri flow 1 x 10 meters is used. 

For late times or when T Tauri is 

5 

not important 5, 10 or 20 x 10 meters 
can be used . 

FINR Fraction of a solar mass being lost 
during the T Tauri stage: 0 ^ FINR. 

We have been using FINR = .5. 

PELT First estimate for a time step (in 

years) in the solution of the thermal 
diffusion equation. PELT must be 
greater than 0 , preferably 10^ years 
for T Tauri heating, 10^ years for 
heating by radionuclides only. 
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MSTEP 


Total number of points in the spatial 
grid for the solution of the thermal 
diffusion equation. It cannot exceed 
200. With DELA = 2 x 10^ meters and 
DELD = 10^ meters MSTEP = 183 and 
IMSTP =83 for the moon. 

IMSTP Number of space points using DELA. 

If r^ is the radial distance from 

the center of the planetary body to 
til 

the l point then 


r 


i 


(i-l)DELA ; 1 s I * IMSTP 

< 

( IMSTP— 1 ) DELA + (i-IMSTP) DELB; 
IMSTP £ MSTEP 


Card Set 5 Name List M4, four entries. 

XI = 0 , XI is the normalized start point for 
the integration of the electric field 
equation. 


INIT = 2 I These two parameters set up the 
NSAV - 0J integration scheme as an initial 

value problem. For cases where the 
internal conductivity of the planet 
becomes too large, INIT is switched 
to zero automatically. NSAV is an 
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internal program parameter which 
must be zero initially but takes on 
the value one for program switching. 
If INIT becomes zero, XI is stepped 
to a value greater than 0 , and the 
field equations are solved as a 
boundary value problem. This is a 
time consuming method and should be 
avoided if at all possible. 

IFct This parameter determines which 

joule heating method will be used. 

1. Unipolar generator (DC limit of 
TM mode ) 

2 . TE mode 

3 . TE mode and unipolar generator 

4 . TM mode 

5 . TE mode and TM mode 

The following data cards are called from the auxilliary 
subroutine CONIN: 

Card Set 6 Name List Cl. 

NT This is the number of times at which 

/ 

thermal profile and output data will 
be printed out. 1 ^ NT ^ 50 
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The time points TB must be set up 
so that TIB (1+1) > TIB (I) 

TIB Time points for output, in years. 

Card Set 7 Name List C2 . 

NT0 Number of points for the starting 

temperature, T0(I) in °C ■ as a 
function of position R0(I) in 
kilometer. 2 ^ NT0 ^ 50 

R0(I) Position at which the temperature is. 

T,0(I) The data from T0 is interpolated 

to cover the points for the thermal 
integration in the main program 
using the subroutine TINTR. 

NSORC This is the number of radio nuclide 
heat sources 1 ^ NSORC ^ 8 

The next card set, set 8, numbering from one to eight, 
depending on NSORL, gives the radio nuclides for the 
planetary heating. The cards are in a 3 x, 3A4, 3E15.5 
format. The 3A4 portion of the format reads in the radio 
nuclide name. 

TEXT (I) Radio nuclide name ( = 1,2,3). 

W(I) Radio nuclide abundance, 

grams/gram. 
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B(I) Radio nuclide heat generation, 

j oules/gram-year . 

FLAM (I) Radio nuclide decay constant, per 
year. 

Card Set 9 Name List C0N. 

DENS The planetary density in grams/cc 

(assumed uniform in this program) . 

FN Index of refraction used in thermal 

conduction formula: FN = 1.7. 

t 

E Energy gap in e.v. for the portion 

of the radiation opacity derived from 
the electronic conduction term: 

E = 3 electron volts. 

\ 

CP Heat capacity: 1.2 Joules/gram-degree . 

FKF Not used in program: = 0 

E0 The opacity (=100/cm) . 

SIG0 Supposed electrical conductivity at 

infinite temperature (=10 mhos/cm) . 

C Lattice thermal conductivity (=7.889 

x 10^ Joule/cm-year- °K) . 
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Card Set 10 Name List C3. 


A1(I) 


1=1, 6. These coefficients are 
for the electrical conductivity 
ct(T) in mhos/meter, where T is 
the temperature. The conductivity 
is given by 


a(T) = ) Al(2i-l)exp <-Al(2i)/T 


M 


In ERAD the conductivity becomes 


a T (T(r)) = a (T) + DP 

where cr(T) is shown above and DP 
is the additive constant in read in 
name list M2 . 

5.5 Test Cases 

A limited number of test cases were run with 
the new computer program. In each case the electrical 
conductivity function used was Olivine (Rikitake 1966) 
with 

o(T) = 0.1e- 5800/T + 10.e- 190 °° /T + 3.5 x 10 V 35000/T 
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Further, only a uniform temperature moon model was used 
with T = 500 °C = 773 °K giving a constant lunar conduc- 
tivity of 5.4 x 10 mhos/meter 

First, consider the TM mode and the unipolar 

generator. For the unipolar generator the joule heating 

-3 

during the T Tauri phase was 8.7 x 10 joules/gram- year 
each from both the radial and tangential electric fields. 
There was a K factor of .37 and a power ratio of 2.4 
(see Sonett et. al . , 1970 for complete explanations of 
these terms) . For the TM mode test case only the n = 1 

_3 

harmonic was used. The resultant heating was 7.4 x 10 

joules/gram- year with a k factor of .36 and a power 

ratio of 2.5. The contributions from the higher harmonic 

_2 

decreases as (2n+l) . The ratio of the TM heating to 

2 

unipolar heating is (8 /t r ) if one takes into account the 
slight difference in the k factors. Note that the 
heating input in this case is four orders greater than 
the radio nuclide heating for a chondritic source,, 

1.28 x 10 - joules/gram-year (Fricker and Reynolds, 1967). 

The contribution from the TE mode is somewhat 
more difficult to analyze. The electric field for the 
TE mode is proportional to the frequency for the, low 
frequencies of interests. Figure 8 shows the TE heat 
input for the first harmonic and the sum of the first 
125 non zero harmonics. Note that the sum term is 
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approximately 100 times greater than the single harmonic. 
This occurs because the amplitude of each harmonic electric 
field is proportional to (2n+l)o> 0 where cd q is the solar 
spin rate and the exciting field decreases as (2n+l) - ^ . 
Thus, except for the affect of attenuation with depth, 
each harmonic produces the same electric heating input. 

The contributions begin to drop off significantly for 
frequencies greater than 0.1 Hertz (27ra/A > 1) . 

5 . 6 Summary 

Only the preliminary portion of the investi- 
gation of planetary heating by time dependent electro- 
magnetic induction has been completed. The major tool 
for the investigation is essentially complete, the computer 
program. Because of the complexity of the program and the 
problem of planetary heating, great care must be taken in 
laying out a plan for examining the possible effects of 
electromagnetic induction on planetary evolution. So far, 
only a constant temperature case has been investigated 
for the TE mode . This case produced an inverted heat 
source with up to ten times greater energy disposition 
than chondritic nuclides toward the planetary surface. 

More conducting moon models would have sharpened the slope. 
Further, in the light of this small data sample, it would 
appear that the TE mode can not be used for the rapid 
heating near time zero. The TM mode heating can deposit 
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10 to 10 more energy. However, the TM mode heating is 
critically dependent on the state of the planetary surface 
electrical conductivity, whereas the TE mode should be 
unaffected by the surface. Therefore, more initial time 
models must be examined before one embarks on a long 
calculation scheme . 
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6.0 


Summary of Report 


6 . 1 Electromagnetic Induction in the Moon 

The major effort for the past year, upon which 
all the other work was based, was the formulation of the 
theory of electromagnetic scattering by a radially inho- 
mogeneous moon. This formulation was a synthesis of 
techniques described by Wyatt (1962) and also used 
originally by Lahiri and Price (1939) and many others. 

The formalism was combined with the lunar boundary con- 
dition first proposed by Blank and Sill (1969) for a two 
layer moon model. The major advance obtained during this 
program was generalization to models where the lunar 
electrical conductivity is an arbitrary function of the 
radius. Discontinuities can also be handled within this 
framework. A spherical harmonic expansion was carried 
out for this theory and harmonics higher than the first 
are retained. The resultant radial wave function has 
been transformed in order to produce a mathematical 
problem which is more tractable for computation; i.e., 
an initial value problem. For special cases, as in the 
near surface thermal gradient calculations of Section 3, 
the solutions were obtained on the computer by a two 
point boundary value calculation which requires much 
longer running time. 
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6.2 


Near Surface Thermal Gradient 


The initial application of the formalism of 
Section 2 was the calculation of the TE transfer func- 
tion for various moon models which could be characterized 
by a constant thermal gradient near the surface. The 
compositional models were derived from a lunar basalt 
(Nagata et. al., 1970) and an olivine (England et. al., 
1968) . The former constitutes one of the most highly 
conducting rock materials while the latter is, for low 
temperature, one of the least conducting rocks. To 
obtain solutions of the differential equations of Section 
2, an iterative technique was carried out to solve a 
boundary value problem. At some inner point, the field 
was assumed zero and the slope of the radial term was 
varied in order to match the surface boundary condition. 

In combination with the initial LSM data (Sonett et. al . , 
1971 a and b) these calculations bracket the thermal 
gradient between 1.4 and 4°K/km. These results were 
based only on the high frequency response, .02 to .04 Hz. 

6.3 Whole Body Lunar Conductivity . 

By using the entire experimental spectrum from 
LSM an iterative fit was obtained for the lunar electrical 
conductivity. The figure of merit for this fit was the 
sum of the squares of the deviation of the calculated 
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transfer function from the experimental transfer function. 
For an eight frequency fit, this value was .06 after five 
iterations. Much more work remains to be done on this 
problem to improve the spatial resolution of the conduc- 
tivity profile. However, a major discovery is the 
apparent presence of a high peak in the electrical con- 
ductivity in the vicinity of r = 1500 kilometers in the 
moon . 

6 . 4 Primordial Planetary Joule Heating 

The mathematical results of Section 2 have been 
incorporated into a computer program for calculating the 
thermal evolution of planetary bodies. The thermal 
sources now include the effects of both the TE mode and 
TM mode. All of the pertinent subroutines have been 
checked for errors and one set of test runs has been 
made for each type of joule heating individually. These 
initial results indicate that the TE mode contribution 
may be significantly greater than the proposed chondritic 
radionuclide sources over the outer lunar shell. Further 
investigation is required, however, to determine 'the 
effects of different thermal and conductivity models on 
this form of heating. The major advantage of TE heating 
is its insensitivity to a low (or zero) conductivity 
surface layer. 
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Appendix A 


The fields (E m , H 111 ) and (E e , H S ) are related 


• rn 0 

to the potentials 0 , 0 according to 


e 1 / a' 


E = E = , 
r r k 


1 ar" 


1 dk l _a_ 2 

k, dr ar + k 
1 / 
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Appendix B 

The spherical harmonic expansions for the 
potential and the electromagnetic fields of the transverse 
electric part of the solar wind excitation are 


Q m = 


M'VH 


sm 





(cos 
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0 ) 
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where P* = t (ui) ' 1 " V 71 


The electromagnetic fields and potential for the trans- 
verse electric mode in the lunar interior are 
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The corresponding formulae for the transverse 
magnetic mode driving field are 
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For the interior fields we have 
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FIGURE CAPTIONS 


1. The electrical conductivity of two representative 

minerals vs. temperatures; a moon basalt after 
Nagata et. al . (1970) and olivine after England et. 

al. (1968) . The former corresponds to a highly 
conductive rock while the latter corresponds to a 
low conductivity rock. 

2. Behavior of the real and. imaginary parts of the TE 
transfer function as a function of radius in the 
moon. This data was obtained using the basalt con- 
ductivity (Nagata et. al., 1970) a = 7.9 e ' 

6 — 1 2 5 /kT 

+5.1 x 10 e ' . The temperature was assumed 

to follow the equation T(r) = 243°K + T (a-r) with 

O 

a = 1740 km . Thermal gradient T of 2, 3 and 
4°K/km were used. The frequency is .04 Hertz. 

3. The TE mode transfer function for the basalt and 

olivine conductivity calculated as a function of the 
thermal gradient for frequencies of .02, .03 and .04 

Hertz . The LSM 



is also shown. 

4. Power spectral density determination for simultaneous 
magnetic field observations on the lunar surface 
(Apollo 12) and in the solar wind near the Moon 
(Explorer 35) (prepared by NASA/ARC) . The y com- 
ponent is east and the z component north at the 
Apollo 12 site. Significant amplification is observed, 
increasing with frequency. Length of the time series 
is 2 hours, 6= 10 gauss. 
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5 . 


The rms tangential lunar transfer function, A = 

2 2 1/2 

QO . 5 (A + A z )J and the normal transfer function 

A^ as a function of frequency. The x direction is 
normal to the lunar surface. The solid line is the 
value of amplification calculated from the conductivity 
profile whose corresponding amplification is fitted to 
experimental values at frequencies of 0.83, 1.7, 5, 

12, 17, 22, 25, and 35 millihz. 

6. Lunar bulk electrical conductivity profiles determined 

from the individual transfer functions A , A , and 
— Y z 

A. The prominent rise of a from the surface inwards 

to r = 1500 km is apparent for all three cases as 

well as the subsequent decrease inwards to r = 1400 km 

followed by a more gradual rise. A tentative version 

of a lunar thermal profile is shown as the grey overlay 

with temperatures indicated on the right hand margin. 

This profile is a fit of conductivities to a Nagata 

basalt in the mantle, an England olivine in the core, 

and the known subsurface temperature of -30°C. The 

insert is shown to suggest extreme values of the a's 

using the one standard deviation limits of the A's 

to calculate conductivities (see Sonett et. al. 1971 

a and b) . 

7. The experimental r.m.s. tangential lunar transfer 
function, A , as a function of frequency. The error 
bars are the same as those in Figures 2 and 4. The 
curve labeled "best fit" are the values of the ampli- 
fication for the conductivity profile labelled A in 
Figure 5. The "2 layer" curve was calculated for a 
model with a core of radius 1560 km and a constant 
conductivity of 7.6x10 ^ mhos/m, and a shell of zero 
conductivity. This value of the core conductivity 
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provides a best fit to A for the given core radius. 

The Dyal and Parkin (1971) amplification curve was 

calculated for a model with an inner core of radius 

_2 

1044 km and a conductivity of 10 mhos/meter, a 

middle layer extending to 1653 km with a conductivity 
-4 

of 1.7x10 mhos/m and an outer layer of zero con- 
ductivity. Other amplification curves are based on 
conductivity models proposed by Ness (1969) and Sill 
(1971) . Note that the "2 layer" curve does not dupli- 
cate the "S" nature of the experimental data. 

T Tauri electrical (joule) heating from the TE mode 
for a 500 °C Rikitake-olivine Moon. The curve 

represents the contribution from the first 125 har- 
monic terms . The TE, curve is the contribution from 
the first harmonic. The horizontal line at 1.25x10 
joules/gram-year is for a uniform chrondritic radio- 
nuclide source after Fricker et. al. (1967) . 
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